Vortices in Trapped Superfluid Fermi-gases 
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We consider a superfluid of trapped fermionic atoms and study the single vortex solution in the 
Ginzburg-Landau regime. We define simple analytical estimates for the main characteristics of the 
system, such as the vortex core size, temperature regimes for the existence of a vortex, and the effects 
of rotation and interactions with normal fermions. The parameter dependence of the vortex core 
size (healing length) is found to be essentially different from that of the healing length in metallic 
superconductors or in trapped atomic BEC in the Thomas-Fermi limit. This is an indication of the 
importance of the confining geometry for the properties of fermionic superfluids. 
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Experiments on cooling of trapped gases of Fermionic 
atoms are at the level where the achievement of 

the predicted BCS-transition ||,|| can be anticipated. 
Cooper-paired trapped atoms would allow to study and 
test the BCS-theory in a controlled manner - for instance 
the classic problem of the BCS-BEC crossover when the 
interparticle attraction varies || could be studied using 
the possibility to tune the interatomic scattering length. 
Several methods for observing the existence of a gap in 
the excitation spectrum of the superfluid Fermi-gas have 
been proposed 

Vortices are a macroscopic signature of the superflu- 
idity, and the vortex core size reflects the typical coher- 
ence lengths of the system. Vortices have played a major 
role in experimental and theoretical studies of superfluid 
Helium ^ and lately, of condensates of Bosonic atoms 
JlO[ . We consider vortex solutions for superfluid trapped 
Fermi-gases. Although a full description of the system 
would require a careful investigation of the Bogoliubov- 
deGennes equations Jll],[l2| , we have chosen as a first at- 
tempt to characterize the system using the Ginzburg- 
Landau equation in a trapped geometry ]l3[ ]. This al- 
lows us to define intuitive estimates for the vortex core 
size and for the temperature regimes where a vortex so- 
lution exists. Furthermore, the effects of rotation as 
well as of interactions with normal fermions (position- 
dependent Hartree-fields) can be incorporated into the 
formalism. We estimate their effects on the critical tem- 
perature and the vortex solution, and develope a self- 
consistent method for treating the effect of the Hartree 
fields. We also discuss how our results can be used in 
estimating the possibilities for observing a vortex. 

The Ginzburg-Landau (GL) equation for atoms in two 
different hyperfine states confined in a symmetric trap of 
frequency O is given by |ll| 
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where A = 2pp \a\/(hn) is the interaction parameter (< 1 
for dilute systems), tJ - 1 is the critical temperature for the 
corresponding homogeneous system, A(i? = t/Rtf) is 
the spatially dependent order parameter and Rtf is the 
cloud size Rtf = (2Ef /(m ti 2 )) 1 / 2 given by the Thomas- 
Fermi (TF) approximation 14|. 

We impose a singly quantized vortex solution for the 
order parameter in cylindrical coordinates in a gen- 
eral anharmonic trap potential: A(R) — A(p,z,8) = 
e 10 A(p, z) where p = , z = and the equation for 
the real order parameter becomes (in h = 1 units) 

1 , d 2 1 d d 2 I 

7 



2m* ^ dp 2 p dp^ dz 2 " 2 



-(p 2 +Az 2 ) 



-E^T)}A(p, z) + C(T)\A(p, z)\ 2 A(p, z) = 0, (2) 



where 



1 48tt 2 (k B T c 

2 7C(3) ^ fl ho 



C 



1 + 2A 

Am* 

7C(3) (Tc 
T 



E„{T) = h 



8tt 2 



(0) 



T 



1 Tc 
= In — 

T 



2^' 



(3) 
(4) 
(5) 
(6) 



and 



3,~. 



T (o) 



ln-^=—. We use the scaling A(R) = 
The asymmetry of the trapping poten- 
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Thus the GL equation becomes analogous to the Gross- 
Pitaevskii (GP) equation for bosonic atoms |l(| for par- 
ticles of an effective mass m* in a harmonic confining po- 
tential with an effective frequency lo* and lj* z = to* VI. 
The nonlinear term plays the role of repulsive interparti- 
cle interaction, and E^(T) is the "chemical potential" of 
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the system. Strong superfluidity, compared to trap en- 
ergy scales, means large effective mass and small effective 
trapping frequency. 

We consider a two-dimensional geometry where the 
single vortex is uniform along its axis (z-axis). Eq.(|^) 
then reads 
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and now E^(T) = \n(T c /T)+ili. We have used the steep- 
est descent method described e.g. in [ fl5| using imaginary 
time propagation to search for the ground state solution 
of Eq.(^). We applied the Crank-Nicholson finite dif- 
ferencing method used to solve the GP equation, but 
without normalization after each step as the size of A is 
determined by the nonlinear term. 

Vortex core sizes The size of the vortex core reflects 
the healing length of a superfluid because within this 
distance, the order parameter "heals" from zero up to its 
bulk value. The first guess to estimate the healing length 
in our case would be to equate the "kinetic energy" term 
in Eq.(0), ~ l/(2m*£ 2 ) to the "interaction energy", as 
was done in case of trapped atomic BEC (GP equation) 
|jl6||. This yields (£ is in Rtf units) 



e = l/(2m*C(T)|A| 2 ), 
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where |A| 2 is the "density" of Cooper-pairs (in ksT c 
units). In the case of the TF approximation for BEC, 
the order parameter is substituded by its value in the 
middle of the trap. Making a corresponding substitu- 



tion here gives £ 



which is the same as the 
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definition of correlation length in metallic superconduc- 
tors in the GL regime. However, it does not correctly 
describe the numerically obtained values for the healing 
length because, for the experimentally feasible parame- 
ters used here, the energy scales in the GL equation do 
not correspond to the TF limit of BEC. 

We propose a measure for the vortex core size by de- 
manding that the kinetic plus potential energy term in 
Eq.(^) has its minimum value. This term corresponds 
to a vortex in a harmonic trap, that is, the first excited 
state. The energetically favoured position is thus the 
maximum of the first excited state wavefunction which 
coincides with the oscillator ground state length. Wc 
thus define 
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Figure I shows the order parameter for a selected set of 
parameters, together with the two estimates Eq.(||) and 
Eq.(||) for the healing length. The results fit excellently 
with Eq.(||) whereas the deviation from Eq.(^) is con- 
siderable and qualitatively different for small and large 



A. The definition of Eq. (j9|) is a function of tempera- 
ture decreasing as T approaches T c (we have confirmed 
this temperature dependence also numerically). This is 
an opposite behaviour to metallic superconductors where 
£ ~ 1/ •v/l — T/T c . The trapping energy becomes rela- 
tively stronger as T — + T c because A decreases. This 
means stronger confinement for A(p) and decreasing £. 
Note that also in BEC, the confinement determines the 
healing length if one is away from the TF regime Jl7| . 
Differences to BEC arise, however, from the temperature 
dependence of the GL equation and from the normaliza- 
tion. 

Critical temperature for a trapped system As pointed 
out in 1 13 1 one can use the GL equation to estimate the 
critical temperature of the trapped system, as compared 
to the corresponding homogeneous case. Close to T c the 
non-linear term is negligible and the GL equation re- 
duces to the Schrodinger equation for a trapped particle, 
E^{T) now denoting the energy. The smallest possible 
energy E^(T C ) = |u) is then simply the ground state en- 
ergy of the trap, lu* + \uu*. Equating |w = lu* + \lu* z 

gives T c in terms of Tc°\ lu* and lu* by the definition 
Cu = |ln(Te /T c ). We use this to calculate the new 
critical temperatures when the effects of rotation and 
Hartree-fields are added. Below we estimate the tem- 
perature regimes for vortex solutions using a similar ar- 
gument. 

Temperature regimes for the existence of a vortex solu- 
tion The ground state energy of the trapped superfluid in 
the quasi-linear regime can be estimated as above to be 
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(in 2-D). A vortex has a higher energy, 



and E^T) = ln(T c /T) + lu* must have large enough 
value in order the GL equation to have a solution. The 
minimum extra energy that a vortex in a trapped non- 
interacting superfluid requires is the 2-D harmonic oscil- 
lator energy lu* (this is also the vortex energy for a non- 
interacting BEC [p^[). Thus the maximum temperature 
T v at which E^(T) can provide this extra energy is given 
by ln(T c /T v ) — lu*. We have checked the validity of this 
estimate for a system of N = 3 * 10 5 atoms, scattering 
lentgh \a\ = 1140A, and trap frequency Q = 820Hz. The 
maximum T at which a vortex solution can exist, esti- 



mated by T c = T v e u = T v e 



0.01464^ 



is T v ~ 0.98T C . The 



actual maximum T where the GL equation gives a nu- 
merical solution was about 0.97T C . The deviation is due 
to neglecting the non-linear term, and T v can be under- 
stood as the upper bound for the maximum temperature. 

Hartree-fields and rotation In practical systems, the 
Cooper-paired atoms are always interacting with the nor- 
mal part of the gas - whose density distribution is now 
position dependent. Moreover, the whole system may be 
rotating, c.f. vortices in atomic BEC |To|]. The GL equa- 
tion used in this paper was derived in Jl3| using the TF 
approximation p4| for the density profile of the trapped 
Fermi-gas in a harmonic symmetric potential. We fol- 
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low this derivation and add a new potential term V(r) 
in the local density approximation. This potential can 
describe e.g. the Hartree fields or rotation. It is assumed 
to be small enough in the sense that one can still use the 
quasi-classical expression for the product of two Green's 
functions in Eq.(7) of [O]. We obtain the GL solution 
as given by Eq.(9) of ]13|] but with a new term in the ex- 
pansion of the Fermi energy around R=0 (see also |20|) 
given by Sf{R)/ef = 1 - R 2 — V(R), where R = v/R T f 
and the Rtf = wp/fi. Note that Rtf and f2 might have 
changed because of adding V(r). Expanding to second 
order and derivating with respect to A* the "potential 
term" in the GL equation now reads 
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where V^- 2 ' (R) denotes an expansion of V to second or- 
der. The potential V(R) has to be smooth enough for a 
second order expansion to be sufficient. 

Effect of rotation on the critical temperature and the 
healing length The TF approximation with a rotation 
term V{p) = —luL z - 



^mu> 2 p 2 gives ef{R)/£f 
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- (2:7) in the rotating frame of reference, where now 
Rtf = v F /^{n 2 -tu 2 ) and 7 = (An 2 /{n 2 - to 2 )) 1 / 2 , 
and A is the assymetry in the trapping potential as de- 
fined earlier. This means that both m* and u>* change. 
As the trapping frequency of the atoms is smaller now, 
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geneous system critical temperature T' c w , and using the 
same kind of procedure as above for a non-rotating sys- 
tem, the new healing length £ r and critical temperature 
can be calcutated. In the limit (^) — > 1, the critical 
temperature — > as for BEC plj , and £ r as well as 
ln(Tc (0) /TT) diverge. 

Effect of the normal fermions The TF approximation 
including the non-paired fermions was introduced in |22f| 
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n(R), where n(R) is 



the density distribution of the atoms in one hyperfine 
state. It was shown that the Hartree-field increases the 
critical temperature of the system. This result is also 
given by our GL treatment: Consider again the weakly 
non-linear regime close to T c . We have a Schrodinger 
equation for a spherically symmetric potential but with 
a new term given by the Hartree field. Because of the 
smooth shape of the fermionic distribution one can as- 
sume n(R) ~ n W(R) and useV^(R) - —^^^n{R) 
in Eq.([10|). Thus we have effectivelly a new harmonic 
symmetric potential m*(uj*) 2 R 2 /2 — C'n(R) ~ —E HF + 
m*(uj* HF ) 2 R 2 /2 with s = (lo* hf /u*) 2 > 1. That is, the 
potential is now deeper and with higher "frequency" . Be- 
cause of the new higher Fermi energy Ep F = Ef + Ehf 
the critical temperature in the corresponding homoge- 



neous system, T c 
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give ln(T c HF{a) /T c HF ) = \n{T^ ] /T^s 1 / 2 . This means 
that T C HF actually deviates more from j^-^O) than T c 
from Tc ; the effect of trapping is enhanced because also 
the normal fermions feel the trapping potential. On the 
other hand, r p^ F< ^> ^ jX ^ anc j f or the parameter values 
we have considered the total effect is that the Hartree 
fields increase the critical temperature (T c f rF ~ 2T C ). 

Self- consistent solution of normal and superfluid 
fermions As a non-rigorous but intuitive first guess to- 
wards a self-consistent treatment, we calculate n(r), in- 
stead of the TF approximation, by the BCS theory in 
local density approximation but with A(r) given by the 
GL equation: 
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where \u k (r)\ 2 , \v k (r)\ 2 = 1/2[1 ± £ k /y/% + \A(r)\ 2 ] and 
^■n(r) — (ji— imfi 2 r 2 ), and / is the 



Fermi-Dirac distribution. We solve n(r) from this equa- 
tion, use it in the GL equation, and iterate until sufficient 
convergence is found. 

The results of the self-consistent calculation are pre- 
sented in Fig. 2. The order parameter increases consider- 
ably, mainly because of the increase in the critical tem- 
perature when the Hartree field is added (T was fixed). 
Quasiparticles fill the vortex core: there is about 5% in- 
crease of n(0) compared to It would be interesting 
to compare our simple self-consistent treatment to the 
rigorous description by the BdG equations. 

Observation of a vortex For the parameters used in our 
calculations, the vortex core sizes/healing lengths vary 
between 2-10 /im. This is close to but still above the 
diffraction limit of light. In principle, for instance the 
laser probing method of l| could be extended to the ob- 
servation of a vortex: the applied Raman beams are fo- 
cused so that they intersect either only in the core, or 
only in the superfluid region, and these two choices give 
absorption peaks at different frequencies. 

Also the smallness of the superfluid fraction makes the 
observation a challenge. In our results, the maximum 
value of the order parameter (vortex height) is of the 
same order of magnitude as the temperature. To esti- 
mate how it depends on temperature and other parame- 
ters, we approximate A max ~ A(p ~ £): we insert p = £ 
into the GL equation, neglect the first derivative (maxi- 
mum) and approximate the second derivative by assum- 
ing a parabolic shape of A(p). This gives (for comparison 
with the numerics, see Fig.l) 
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increases. Simple considerations 



Note that our estimate for the maximum vortex temper- 
ature T v corresponds to A max being real. 

We have used the Ginzburg-Landau equation in a 
trapped geometry to define analytical estimates for the 
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basic quantities describing a trapped superfluid Fermi- 
gas: the vortex core size and height, maximum temper- 
ature for a vortex solution, and the changes caused by 
additional potentials such as rotation or spatially varying 
Hartree fields. A striking difference to metallic supercon- 
ductors was found in the temperature and system param- 
eter dependence of the vortex core size/healing length. 
Our results indicate that the effect of the confining ge- 
ometry is essential for mesoscopic fermionic superfluids, 
especially when considering excited state solutions such 
as vortices. 
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FIG. 1. Vortex solutions for the order parameter in 
fcflTc units and in trap units. Solid lines are for 
n = 2Tvl00Hz,lkHz,3kHz (N = 3 * 10 5 , |a| = 1140A, 
T = 0.89T C ). Dashed lines correspond to N = 10 5 , 3* 10 5 , 10 6 
(T = 0.8T C , n = 820Hz and \a\ = 1140A). Dot-dashed lines 
are for |o| = 985, 1118, 1608A (TV = 3 * 10 5 , Vl = 820#z, 
T — 0.89T C ). In all the three cases, the highest value of Q, 
N or | a. | corresponds to the curve with largest maximum A. 
The healing lengths given by Eq.(^) and Eq.(^|) are repre- 
sented by * and +, respectively. The estimate for A max given 
by Eq.(|l2|) is represented by o. 




FIG. 2. Vortex solutions for the order parameter are given 
in ksT c units and density of normal fermions in one hy- 
perfine state in 10 12 atoms /cm 3 . Solid lines are the inde- 
pendent initial solutions for Eq.(ll) and the GL equation, 
the dot-dashed lines are the self-consistent solutions. Here 
T = 0.9T C , N ~ 3 * 10 5 , n = 820Hz and \a\ = 1140A. 
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